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This  thesis  identifies  the  mean  orbital  elements  which  produce  arcs 
of  minimum  altitude  variation  over  an  oblate  planet  with  an 
axi- symmetric  gravitational  field.  Such  orbits  are  useful  for 
surveillance  or  scientific  study  missions  using  optics  with  fixed  focal 
lengths . 

Both  Earth  and  Mars  are  considered  and  the  optimum  eccentricity  is 
found  as  a  function  of  argiiment  of  periapsis  and  inclination  for  two 
values  of  semi -major  axis  for  each  planet.  The  results  are  curve  fit 
to  develop  a  single  equation  which  identifies  the  eccentricity  needed  to 
produce  an  arc  of  minimum  altitude  variation  given  the  argument  of 
periapsis,  incliiiation,  semi-major  axis,  ellipticity  of  the  planet, 
equatorial  radius,  and  the  zonal  J^. 

Once  arcs  with  minimum  altitude  variations  are  identified,  the 
properties  of  the  arcs  are  considered.  The  mid-latitude,  altitude, 
duration,  and  latitude  range  of  the  arcs  are  found  as  a  function  of 
argument  of  periapsis  and  inclination  for  various  planet  and  semi-major 
axis  combinations. 

The  secular  change  in  mean  orbital  elements  is  considered  to 
determine  the  most  stable  orbits.  Secular  changes  in  orbital  elements 
due  to  the  geopotential,  drag,  and  third  body  effects  are  considered. 

The  velocity  impulse  needed  to  return  the  satellite  to  the  original 
orbit  from  the  perturbed  orbit  is  found  and  used  to  determine  stability. 
Identifying  orbits  which  require  minimum  station  keeping  fuel  allows 
planners  to  select  orbits  permitting  longer  useful  operational  life. 
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ORBITS  CONTAINING  ARCS  OF  MINIMUM  ALTITUDE  VARIATION 


I .  Introduction 


Background 

Mean  orbital  elements  may  be  determined  to  produce  orbits  about  an 
oblate  planet  which  contains  arcs  of  constant  altitude[6].  The  arcs  of 
constant  altitude  may  be  positioned  over  specified  latitude  ranges  -- 
this  ability  is  useful  for  surveillance  or  scientific  study  missions 
using  optics  with  fixed  focal  lengths.  This  thesis  focuses  on  finding 
orbits  with  arcs  of  constant  ellipsoidal  altitude  for  axi-symmetric 
planets  (only  zonal  harmonics  are  considered  in  the  geopotential). 
Because  of  these  restrictions,  the  longitude  of  the  ascending  node  is 
not  considered  as  a  variable.  Given  the  required  orbit,  arcs  of  nearly 
constant  altitude  may  be  identified.  The  eccentricity  becomes  a 
function  of  the  required  inclination,  argument  of  periapsis,  semi-major 
axis,  and  planet  parameters  (ellipticity ,  equatorial  radius,  and 
geopotential).  The  independent  variables  are  also  functions  of  the 
desired  altitude  and  latitude  of  the  arc. 

Mean  orbital  elements  which  produce  stable  orbits  need  to  be 
considered.  Missions  requiring  these  types  of  orbits  will  typically 
have  low  altitudes,  so  drag  needs  to  be  addressed.  Kalil  [6]  presents 
a  first  order  analysis  neglecting  drag.  Drag  primarily  affects 
eccentricity  and  semi -major  axis;  consequently,  it  can't  be  neglected 
for  missions  longer  than  one  orbital  period.  The  change  in  semi -major 
axis  due  to  drag  will  act  to  reduce  the  altitude  of  the  arc  of  constant 
altitude  in  a  few  orbits.  Eventually,  reducing  the  semi-major  axis  will 
result  in  an  orbit  which  no  longer  contains  an  arc  of  constant  altitude; 
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however,  the  orbit  is  much  more  sensitive  to  a  change  in  eccentricity, 
inclination,  or  argument  of  periapsis  --as  the  drag  reduces  the 
eccentricity,  the  resulting  orbit  will  no  longer  contain  an 
arc  of  constant  altitude. 

Definition  of  Constant  Altitude  Arc 

A  constant  altitude  arc  is  defined  for  the  purpose  of  this  thesis 
as  a  portion  of  an  orbit  during  which  the  ellipsoidal  altitude  of  the 
satellite  changes  only  slightly  with  time. 


Ob  ject ive 

Mean  orbital  elements  of  orbits  containing  constant  altitude  arcs 
are  identified  for  both  Earth  and  Mars.  The  planets'  ellipticity,  zonal 
harmonics  up  to  the  third  order,  atmospheric  drag,  and  third  body 
effects  are  considered.  For  a  given  planet,  the  required  orbital 
elements  are  identified  to  achieve  a  desired  altitude  at  a  given 
latitude.  Secular  changes  in  the  mean  orbital  elements  are  then 
considered  in  order  to  locate  the  most  stable  orbits. 

Approach 

A  search  method  is  employed  to  determine  the  eccentricity,  as  a 
function  of  argument  of  periapsis  and  inclination  for  several  cases 
(varying  the  following  parameters;  planet,  semi -major  axis,  and  the 
perturbative  accelerations  considered).  First,  the  two-body  problem  is 
considered  and  an  empirical  relationship  is  developed  which  determines 
eccentricity  as  a  function  of  argument  of  periapsis  and  a 
non-dimensional  grouping  of  ellipticity,  equatorial  radius,  inclination,' 
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and  semi -major  axis.  Next,  additional  perturbative  accelerations  are 
addressed  [1,2,8.13,17],  For  the  short  term  considered  (one  period), 
the  J2  term  of  the  geopotential  causes  the  most  significant  change.  An 
empirical  relationship  to  predict  the  change  in  eccentricity  from  the 
two-body  results  is  determined  as  a  function  of  argument  of  periapsis 
and  a  non-dimensional  grouping  of  J2,  inclination,  semi-major  axis,  and 
the  equatorial  radius. 

The  choice  of  w,  I,  and  a  which  minimizes  secular  variations  is 
then  considered.  Secular  perturbations  due  to  third  bodies]?],  zonal 
harmonics  up  to  sixth  order[12],  and  atmospheric  drag  (including  a 
rotating  atmosphere [ 16 ] )  are  addressed.  The  sum  of  the  velocity 
impulses  needed  to  correct  the  secular  changes  in  the  orbit  caused  by 
perturbations  is  calculated  and  used  in  identifying  the  most  stable 
orbits . 
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II.  Theory 


In  order  to  determine  orbital  elements  of  orbits  which  contain  arcs 
of  minimum  altitude  variation,  the  motion  of  the  satellite  must  be  know. 
Kepler  identified  satellite  position  and  velocity  under  the  influence  of 
a  central  body  modeled  as  a  point  mass.  This  two-body  solution 
identifies  the  position  of  the  satellite  using  the  six  classical  orbital 
elements  (semi-major  axis  (a),  eccentricity  (e),  inclination  (i), 
argument  of  periapsis  (w) ,  longitude  of  the  ascending  node  (0),  and  true 
anomaly  (v))  .  Perturbations  from  two-body  motion  due  to  small 
accelerations  may  be  calculated.  Averaging  the  change  in  orbital 
elements  over  one  period  eliminates  periodic  short  term  variations  [15]. 
The  secular  changes  in  orbital  elements  identify  changes  which  must  be 
corrected  to  return  the  satellite  to  the  original  orbit.  This  section 
will  present  two-body  motion  and  perturbations  from  the  two-body  motion 
due  to  non- spherical  gravitational  fields,  third  bodies,  and  atmospheric 
drag.  The  secular  changes  in  classical  orbital  elements  will  be 
presented  for  each  perturbing  acceleration  considered. 

The  equations  of  motion  presented  in  this  section  will  be  evaluated 
with  various  combinations  of  classical  orbital  elements  to  Identify 
combinations  of  elements  which  produce  orbits  with  minimum  altitude 
variation  arcs. 

The  secular  variations  of  classical  orbital  elements  presented  in 
this  section  will  be  evaluated  for  orbits  containing  arcs  of  minimum 
altitude  variation  and  will  be  used  to  determine  the  velocity  impulse 
required  to  maintain  the  satellite  in  the  desired  orbit. 
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Equation  (1)  may  be  solved  to  determine  the  position  and  velocity 
of  a  satellite 

•  • 

r  =  -W  =  a  (1) 

where  a  is  the  sum  of  accelerations  acting  on  the  satellite  and  V  is  the 

potential  function.  The  two-body  solution  and  perturbations  from 

two-body  motion  due  to  a  non-spherical  gravitational  field,  third  body, 

and  drag  will  now  be  presented. 

Two-Body  Motion.  For  two-body  motion,  equation  (1)  may  be 

2 

rewritten  by  recognizing  W  =  p/r 

r  +  ^  -  0  (2) 

r 

where  r  is  the  position  of  the  satellite,  and  //  is  the  product  of  the 
universal  gravitational  constant  and  the  mass  of  the  planet.  Equation 
(2)  may  be  solved  to  provide  the  satellite's  distance  from  the  center  of 
the  planet.  Equation  (3)  presents  Kepler's  equation: 

r  «=  a(l  -  e^)/[l  +  e  cos(v)]  (3) 

Where  r  is  the  distance  from  the  center  of  the  planet,  a  is  the 
semi -major  axis  of  the  orbit,  e  is  the  eccentricity,  and  v  is  the  true 
anomaly . 

Geopotential .  The  geopotential  function  for  a  spherical 
gravitational  field  is  -p/t  The  geopotential  function  for  a 
non-spherical  gravitational  field  is  presented  in  equation  (4)  [1] 


00  00 


'"El 


£m 


i*0  m*0 


j  a 


im 


-f+i  2_  I  _  „  , 

p-0 


sinff  i-m  even 
sin//  i-m  odd 


(4) 


(5) 


where 


tmax 


=  I  - 

t!  (i-t)! 


(2i-2t)! 


Ill 

I  (" 

s=0 


s . 
cos  i 


(i-m-2t)! 

gu 

g-k  fi-m-2t+s 


.  i-m-lt  . 
sin  1 


XcD-''  ( 

g=Si 


g 


]  (p";-J 


(6) 


=  max[0,p-t-m-s] 
g^  =  min[  i-ni-2t+s,p-t] 

H  =  (l-2p)(ca+j/)  +  m(n+S  ) 

O 


'  p  ,  p<k 

tmax  -  * 

.  k  ,  p>k 


k  - 


By  considering  only  zonal  harmonics,  equation  (6)  may  be  simplified: 


V  - 


iO 


(7) 
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f^K 


iO 


i+l  "'iO 


i-2t 


I  - 

t!  (i-t)! 


(2i-2t)! 


(i-2t)!  2 


2(i-t) 


.  i-2t. 
sin  1 


Y  f  ,  k+c  r  cos(2-2t-2c)  (i4>+u)'\  2  even 

^  [  c  J  ^  sin('i-2t-2c;('w+i/;J  i  odd 


c=0 


where 


Recognizing 


cosx  =  cos-x 
-sinx  =  sin-x 

and  defining  as  -C^q  results  in  further  simplification  of 
equation  (8). 

J  k 


mr: 


J 

i+l  2 

o 

rt 

1 

ft 

r 

’  t] 

(-1)’' 

X 

2  pp+l] 
L  1  P  J 

+  2  )  pP 

(-1)^  sin(u>+v) 

.  i-2t. 
sin  I 


(-1) 


k+c 


cos[  (2p-2c  )(ijf+i/)] 

sin/'  (2p-2c+l)  (u)+v) ] 


(9) 


where  the  top  choices  in  equation  (9)  are  used  when  2  is  even,  bottom 
choices  are  used  when  2  is  odd,  and  p  -  k-t 

Merson  [12]  has  averaged  the  perturbations  due  to  zonal  harmonics 
over  one  period,  his  results  are  presented  in  equations  (10)  and  (11). 
The  equations  for  change  in  argument  of  periapsis  (w)  and  longitude  of 
ascending  node  (0)  are  not  valid  if  the  inclination  (i)  is  zero;  also, 
the  equation  for  change  in  argument  of  periapsis  is  not  valid  if  the 
eccentricity  (e)  is  zero. 
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For  semi-major  axis: 


Aa  -  -  (18/4)jrJ2^R^a  ^(1-e^)  ^esinu>(l+ecosw)^(4-5sin^i)  (10) 


Changes  in  the  other  orbital  elements  are  described  by  equation  (11) 
where  Ax  represents  either  Ae ,  Ai ,  Aw,  or  AfJ;  x^  represents  e^^,  ij^,  Wj^, 

or  (k  =  2, 3, 4, 5, 6);  and  X22  represents  £22.  ^22’  ‘^22'  ^22' 

example,  to  calculate  the  secular  change  in  eccentricity  over  one 
period,  replace  Ax  with  Ae ,  Xj^  with  ej^,  and  X22  with  622  equation 
(11).  The  equations  used  to  calculate  the  variables  in  equation  (11) 
are  presented  below  for  each  orbital  element. 


Ax  -  2^[XJj^(R/P)*"xj^  +  J2^(R/p)%2  5 


(11) 


For  eccentricity: 


®2 

®3 

®5 

®6 

®22 


=  0 

3  2  2 

=  -  g  (1-e  )  sini  cosw  (4-5sin  i) 

=  -  ^  (l-e^)(6-7sin^i)sin^i  e  sin2w 

=  2^  (l-e^)sini[ (8-24sin^i+21sin^i)(8+6e^)cosw 
+  7(8-9sin^i)sin^i  e^  cos3w] 

-  (l-e^)sin^i[ (16-48sin^i+33sin^i)(10+5e^)sin2w 
+  3(10-llsin^i )sin^i  e^  sin4w] 

-  sinw[6(4-5sin^i)(l+ecosw)^ 


-  sin^i(l-e^) { (16-20sin^i)-(14-15sin^i) }ecosw] 
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For  inclination: 


^3  ”  "  8  ^  cosi  cosw  (4-5sin  i) 

,  45  2,,  .  2.,  .  2.  .  „ 

I4  =  -  ®  (6-7sin  i)sin  x  sin2t.j 

i^  =  -  e  cosi [ (8-24sin^i+21sin^i) (8+6e^)cosw 

2  2  2 

+  7(8-9sin  i)sin  1  e  cosSoj] 
is  =  -  J'^2Q  ®  sin2i  [(16-48sin  i+33sin  i)(10+5e‘^)sin2w 

+  3(10-llsin^i)sin^i  sin4w] 

3  2 

^22  *  6^  sin2i  [(32-40sin  i)e  sinw) 

+  ( -14+15sin^i)e^  sin2w] 

For  longitude  of  the  ascending  node  (i»<0): 

^2  ”  -(3/2)cosi 

0^  ”  (3/8)  (4-15  sin^i)  e  sino;  coti 

^  cosi [ (4-7sin^i)(6+9e^)-3(6-14sin^i)e^cos2w] 

1  c  0/0 

^5  ”  256  i+33sin  i)(8+6e^)e  sinw 

+  7(8-15sin^i)e^  sin3w] 

Og  -  -  ^0^240  [  10(8'36sin^i->'33sin^i)(8+40e^+15e^) 

-  25(16-96sin^i+99sin^i ) (4+2e^)e^cos2w 

-  15(20-33sin^i)sin^i  e^  cos4<.>] 

^22  “  (3/32)cosi[ (12-80sin^i)  +  16(4-10sin^i)ecosw 
+  ( -4-5sin^i )e^  +  2( -7+15sin^i)e^cos2w) 
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For  the  argument  of  periapsis  (i»‘0,  e^'O): 


<*>2  -  (3/4)(4-5sin^i) 

-1  0  7  2  2 
-  (3/8)e  sinwsini  [  (4-5sin  i)  +  (35cos  i  -  4cosec  i)e  ] 

(16-62sin^i+49sin^i)  +  (6sin^i-7sin^i)cos2w 

^  /1Q  •  2.  ^  189  .  4.,  2 

+  (18  -  63sin  1  +  -2|— sin  i)e 

+  (-6  +  35sin^i  -  ^sin^i)e^cos2u>] 


re  sinwcoseci 


2  3  4 

+  2sin  i  -  jsin  i 


1  •  2. 
sin  1 


9.2. 


+  -1  +  gSin  ilsin  i  e  cos2o> 


105  -1  .  rf  ^  ^  o  •  2.  3  .  4.1  . 

=  “16®  sinwcoseci  I  -  j  +  2sin  i  -  jsin  ilsii 

^  f4  87  .  2.  ^  67  .  4.  357  .  6.1  2 

+  jy  -  —Sin  1  +  — ysin  1  -  — jgSin  ije 

f  9  2  1  4  2 

+  |*1  +  8®^^  ilsin  i  e  cos2u 

,  (I  .  ^  267^.„«.  . 

+  l^l  -  ^sin^i  +  ^sin^ijsin^i  cos2u>j 
525  r8f,  _  .  2.  129  .  4.  297^.  6.1 

‘"6  “  "  ■  32®""  "J 

(0  /  0 
2  -  6sin^i  +  — gsin^i  sin'll  cos2u>j 

+  6(l  - 

+  ^-2  +  25sin^i  -  ^sin^i  +  ^sin^ije^cos2w 

^  3f  ,  11  .  2.1  ,  4.  2  . 

+  gl  1  -  Yo®""  -ilsin  X  e  cos4w 

*  2(-  ^ 

+  [-1  +  ^sin^i  -  ^||sin^i  +  ^||sin^ije^cos2w 
+  +  ^sin^i  -  ^^sin^ijsin^i  cos4wj 


fO  /  0  1 

2  -  6sin^i  +  — gsin^i  sin'll  cos2u>J 
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2.  5.4. 

sin  1  -  -^sin  1 


COS(j) 


f95  .  2.  445  .  4.1  f  o  ^  23 

1  -  -^sin  -iJ  +  [-2  +  Y2 


.  2.  ^  5  .  4.',  „ 

sin  1  +  gSin  i|cos2w 


O' 


r  25  ^  461  .  2.  50  .  4.1  1  ^  5  .  2.1  , 

I — ^  +  — ^sin  1  -  —sin  ijecosw  2  ■*■  ijecosjw 


f  7  3.2.  15  .  4.1  2^f  7  79  . 

[12  ■  8®^"  ^  ■  32®^"  '^[12  '  24®^" 


2.  ^  45  .  4.1  2  , 

1  +  YgSin  lie  cos2lo 


J 


These  equations  will  be  used  to  determine  the  secular  changes  in 
orbital  elements  due  to  perturbations  caused  by  the  geopotential.  The 
changes  will  be  combined  with  secular  changes  due  to  third  bodies  and 
atmospheric  drag  to  determine  the  velocity  impulse  required  to  keep  a 
satellite  in  the  desired  orbit  which  contains  a  constant  altitude  arc. 

Third  Body.  For  Martian  cases,  perturbations  due  to  the  Sun  will 
be  evaluated,  while  perturbations  due  to  both  the  Sun  and  Moon  will  be 
evaluated  for  Earth  orbits.  The  analysis  presented  below  is  applicable 
to  any  third  body.  Figure  1  demonstrates  the  relationship  between  the 
central  body,  the  satellite,  and  the  third  body. 


Third  Bodg 


Figure  1.  Geometry  for  Third-Body  Perturbations 

The  position  of  the  satellite  with  respect  to  the  central  body  is  r,  r' 
is  the  position  of  the  third  body  with  respect  to  the  central  body,  and 
p  is  the  position  of  the  third  body  with  respect  to  the  satellite.  The 
angle  between  r'  and  r  is  S,  and  p  -  r'  -  r 
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Equation  (12)  presents  the  disturbing  function  due  to  a  third  body 
where  the  perturbing  accelerations  due  to  the  third  body  are  VR' . 


R' 


(12) 


Recognizing 


and 


P  =  (r'-r)-^ 

=  r'[l  +  (r/r*)^  -  2(r/r ' )cosS] 

r ' • r  =  r ' rcosS 


results  in  equation  (13),  an  expansion  of  equation  (12) 


R' 


rcosS  1 

r’  J 


(13) 


The  binomial  expansion  is  employed  to  expand  the  first  term  in 
equation  (13) 


[1  +  (r/r')^ 


2(r/r')cosS]‘-^  ["„• ^] ( -1)"[ (r/r ' )^- 

n=0 


2(r/r '  )cosS]’^ 


Neglecting  terms  of  (r/r')  equal  to  or  greater  than  third  order  results 
in  an  approximation  to  equation  (12).  For  the  case  of  Earth's  moon,  if 

3 

the  satellite  orbit  remains  within  six  Earth  radii,  (r/r')  <  .001.  For 

perturbations  due  to  the  Sun  this  term  will  be  much  less  than  unity  for 
both  Earth  and  Mars.  Equation  (14)  incorporates  these  approximations  to 
equation  (12) 
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2 

R'  =  rr  +  (3cos^S  -  1)  (14) 

r  2r'^ 

2  3 

Substituting  =  n'  a'  in  equation  (14)  and  observing  V(/i  3/r')-0 
results  in  a  simplified  potential  equation  for  a  third  body 

R'  =  (a^n' V2)(r/a)^(a'A' )^(3cos^S  -  1)  (15) 

Kaufman  [7]  has  developed  the  secular  change  in  orbital  elements  due  to 
third  body  perturbations.  Equations  (16)  -  (20)  are  the  averaged  change 
in  orbital  elements  for  one  period.  Equations  (18)  and  (19)  are  not 
applicable  when  the  inclination  is  zero. 


Aa  =  0 

Ae  =  (-15jr)  (l-e^)'^  (n'/n)^  ea/8(aVr')^ 

Ai  =  -6ir  (n'/ti)^  (1-e^)  (sini)  ^  (a'A')^ 
X  [aS(l+4e^)  +  j0«(l-e^)  -  Sa^e^cosi] 

An  =  6jr  (n'/n)^  (l-e^)"'^  (sini)’^  (a'/r’)^ 

2  2 

X  [a7(l+4e  )sinw  +  ^7(l-e  )cosw] 


(16) 

(17) 

(18) 

(19) 


Aw  =  67r  (n'/n)^  (a'/r')^  (l-e^)'^  ^20) 

X  {(4a^-;0^-l)  -  (7COsi/sini)(l-e^) [Q(l+4e^)sinw  +  ^(l-e^)cosw] ) 


where : 


P 


Q 


_  cosO'cosit' 

r'°  -  sinfl'cos^' 
sinfl'cos^'' 

cosOcosw  -  sinflsinwcosi 
sinflcosw  +  cosOsinwcosi 
sinicosw 


-cosOcosw  -  sinfisinwcosi 
-sinflcosu  cosHsinwcosi 
sinicosw 


_  sinnsini 

R  -  -cosOsini 
cosi 
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P' 


-sinflcosw  -  cosOsinwcosi 
cosOcosw  -  sinflsinwcosi 
0 


Q' 


sinOcosu)  -  cosfisinwcosi 
-cosOcosw  -  sinnsinwcosi 
0 


a 


r'° 


•P 


•Q 


7 


=  r'°R 


S 


7'° 


•P' 


£ 


•Q' 


Equations  (16)  through  (20)  will  be  used  to  calculate  the  average 
change  in  orbital  elements  over  one  period  due  a  third  body. 

Atmospheric  Drag.  The  acceleration  due  to  drag  acts  in  the 
opposite  direction  from  the  velocity  of  the  satellite  with  respect  to 
the  atmosphere  and  is  modeled  by  equation  (21): 


where 


D  = 


(21) 


/S  -  CpS/2m 


and  is  the  drag  coefficient,  S  is  the  cross-sectional  area  of  the 
satellite  normal  to  the  velocity  of  the  satellite  with  respect  to  the 
atmosphere,  m  is  the  satellite  mass,  p  is  the  atmospheric  density,  and 
Vg  is  the  velocity  of  the  satellite  relative  to  the  velocity  of  the 
atmosphere . 

The  atmospheric  density  at  a  given  location  doesn't  remain 
constant,  but  changes  due  to  the  time  of  day,  solar  activity,  and 
season.  Consequently,  modeling  the  atmosphere  accurately  is  extremely 
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difficult  and  the  best  results  are  obtained  by  employing  models  or 
tabulated  mean  values  which  vary  with  time.  However,  the  atmospheric 
density  may  be  approximated  as  decreasing  exponentially  with  increasing 
altitude  (see  equation  (22)). 


P 


^0^ 


-z/H 


(22) 


where  is  atmospheric  density  at  some  reference  height  above  the 
surface  of  the  planet,  z  is  the  actual  altitude  minus  the  reference 
altitude,  and  H  is  the  scale  height  of  the  atmosphere. 

The  velocity  term  in  equation  (21)  is  the  velocity  of  the  satellite 
with  respect  to  the  atmosphere.  A  technique  to  approximate  the  velocity  in 
terms  of  orbital  parameters  is  presented  [19:275-276].  Let 


V  =  V  +  V 
s  a 

where  V  is  the  velocity  of  the  satellite  with  respect  to  an  inertial 
reference  frame  with  its  origin  at  the  center  of  the  planet,  is  the 
velocity  of  the  satellite  relative  to  the  velocity  of  the  atmosphere, 
and  is  the  velocity  of  the  atmosphere  with  respect  to  an  inertial 
reference  frame  with  its  origin  at  the  center  of  the  planet. 

Assuming  the  angular  velocity  of  the  atmosphere  (w  )  is  uniform 
about  the  north-south  axis  results  in  an  atmospheric  velocity  in  terms 
of  the  geocentric  latitude  and  satellite  radius  (r). 


V  -  rw  cosif> 
d  s 


If  \l>'  is  the  angle  between  V  and  then 


(24) 
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(25) 


V  ^  +  V  ^  -  2W  cosV>' 

s  a  a 

The  atmospheric  velocity  vector  is  in  the  local  horizontal  plane  and  the 
satellite's  velocity  vector  is  almost  in  the  local  horizontal  plane. 

The  angle  between  the  two  vectors  (\l>' )  may  be  approximated  by  the  angle 
between  the  projection  of  the  satellite's  velocity  vector  onto  the  local 
horizontal  and  the  atmospheric  velocity  vector  (V’)[9].  Consequently, 
using  the  law  of  cosines  from  spherical  trigonometry, 

cosV>  cos^  -=  cosi 

V  costj)'  =  V  cosrj>  -=  rw  cosi  (27) 

d  3  d 

The  satellite's  velocity  relative  to  the  atmosphere  may  be  approximated 
in  terms  of  its  velocity  with  respect  to  an  inertial  reference  frame, 
orbital  elements,  and  angular  velocity  of  the  atmosphere. 

V  ^  +  V  ^  -  2W  cosib 

s  a  a 

2  2  2  2 

-  V  +  r  w  cos  ^  -  2Vrw  cosi 
a  a 

=  V^[l  -  (rw  /V)  cosi]^  +  r^w  ^(cos^0  -  cos^i) 

d  d 

2  2  2  2 
V  may  be  further  approximated  by  recognizing  r  w  <  .005  V  [9]; 

V  ^  =  V^[l  -  (rw  /V)cosi]^  (28) 

s  d 

Consequently, 

D  -  p/9V^[l  -  (rw  /V)  cosi]^  (29) 

d 

Sterne  [16]  has  considered  secular  changes  in  the  classical  orbital 
elements  due  to  a  rotating  atmosphere  about  an  oblate  planet.  Equations 
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(30)  -  (34)  present  the  secular  changes  in  orbital  elements  over  one 
period. 


X  [cos£  -  .  5d(Ej^ )  (2cos£-e-ecos^£)/(l-e^ )  ]  dE 


(30) 


(31) 


Ai  =  -  .  5(«^/n)/9asini ( 1 -e^) 


.  ,4^1 


X  {1  +  cos2coE^  ^[(2-e^cos^£  -  1  +  2e^  -  2ecos£]  )  dE 


(32) 


2?r 

Afl  =  -  .5(w  /n)/9asin2w(  1 -e^)  f  p  (l-e^cos^£) 

^  oJ 

X  [2e^  -  1  -  2ecos£  +(2-e^)cos^£]  dE 


(33) 


Aw  =  -cosi  AO 

where 


Ej^  -  (1  -  ecosE) 

^2  “  (^  ■*■  ecosE) 
d  -  (w  /n)cosi(l-e^)^'^^ 

O 


^  /  3.1/2 
n  -  (M/a  ) 


(34) 


Equations  (30)  through  (34)  will  be  numerically  integrated  to  determine 
the  average  change  in  orbital  elements  over  one  period  due  to  drag. 
Secular  changes  in  orbital  elements  due  to  all  the  perturbations 
considered  will  be  combined  to  determine  the  total  change  in  orbital 
elements  over  one  period.  For  orbits  about  Earth  and  Mars  perturbations 
due  to  the  geopotential,  atmospheric  drag,  and  the  Sun  will  be 
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calculated.  For  Earth  orbits,  the  secular  changes  in  mean  orbital 
elements  due  to  the  Moon  will  also  be  calculated. 


Planet  Geometry 

In  this  section,  the  radius  of  an  oblate  planet  is  found  as  a 
function  of  latitude.  The  radius  is  needed  to  calculate  the  ellipsoidal 
altitude  which  is  the  difference  between  the  radius  of  the  orbit  and  the 
radius  of  the  planet. 

Consider  an  axi- symmetric  planet  with  semi-major  axis  (a), 
semi-minor  axis  (b) ,  and  radius  (R)  as  a  function  of  latitude  (^).  The 
ellipticity  (f)  of  the  planet  is  defined  in  terms  of  the  semi-major  (a) 
and  semi -minor  (b)  axes  as 


f  -  1  -  b/a 

or  in  terms  of  the  planet's  eccentricity  (e) 

f  =  1  -  (1  -  e^)-^ 

A  relationship  to  determine  the  planet's  eccentricity  is  presented 
below. 


2  2 
1  -  e^  -  (1  -  f)^ 

Equation  (36)  may  be  employed  to  approximate  the  planet's  radius  given 
the  latitude  and  the  ellipticity: 


R  -  a[l  -  f  sin^^  +  0(f^)] 


(36) 
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The  reasoning  employed  to  achieve  this  result  is  presented.  Consider  an 
ellipse  circumscribed  by  a  circle  (see  Figure  2).  The  geometry  of  the 
ellipse  dictates  the  relationship  between  the  semi -major  and  semi -minor 
axes,  the  eccentricity,  and  any  point  on  the  perimeter  of  the  ellipse. 


u 


Figure  2 .  Geometry  of  Planet 


1  =  (x/a) 
.  2 

ae  =  (a 


‘  +  (y/b) 

b2)'5 


2 


Eventually,  an  expression  for  R/a  as  a  function  of  ^  must  be  found. 
2 

First,  (R/a)  is  found  as  a  function  of  E. 

o2  2^2 

R  =  X  +  y 

“  [a  cosE]^  +  [b  sin£]^ 

-  [a  cosE]^  +  [a  (1  -  e^)'^  sinE]^ 

2  2  2  2  2 

-  a  [cos  E  +  sin  E  -  e  sin  E] 

2  2  2 
(R/a)^  -  1  -  e^sin^E 
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Sin£  is  found  as  a  function  of  ^  and  length  y  may  be  expressed 
in  terms  of  the  semi -major  axis,  eccentricity,  and  sin£: 

y  =  R  sin^ 

■=  b  sinE 

=  a  (1  -  e^)'^  sinE 

These  relationships  are  used  to  solve  for  sinE. 

sinE  =  (R/a)  [  sin<^/(l  -  e^)'^] 


Consequently , 

(R/a)^  =  [1  +  F  sin^(^]'^ 


where 

F  -  eVd  -  e^)  -  l(l-f)'^  -  1] 


Equation  (37)  presents  a  relationship  to  evaluate  the  radius  given 
the  eccentricity  and  latitude. 


R  =  a 


2  2 

1  +  e  cos  4> 


(37) 


Equation  (37)  will  be  approximated  by  a  truncated  series  expansion.  The 
series  will  be  developed  by  expanding  (R/a)  and  F  using  the  binomial 
expansion.  Equation  (38)  presents  the  expansion  of  F. 

F  -  -I  . 

n-O 
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(38) 


00 


F  ^  I  (n+1)  f" 
n-1 


Equation  (39)  presents  the  expansion  of  (R/a); 


R/a  =  tF  sin^^]" 

n=0 


(39) 


Combining  equations  (38)  and  (39)  results  in  equation  (40)  which 
determines  the  radius  of  the  planet  as  a  function  of  latitude, 
ellipticity,  and  equatorial  radius. 


R  =  Re^['^n^]  [  ^(“+1)  f”]""  [sin%]  (40) 

n=0  m=l 

Truncating  third  order  terms  and  higher  from  equation  (40)  results  in 
an  approximation  of  the  radius  as  a  function  of  latitude,  ellipticity, 
and  equatorial  radius : 

R  -  R^  [1  -  f  sin^^  +  1.5  f^  (sin^^  -  sin^4>)  +  0(f^))  (41) 

4  2 

Consider  the  term  (sin  4>  *  sin  4>)  in  equation  (41).  Taking  the 
derivative  with  respect  to  4>  indicates  optimiims  will  occur  at  0°, 
±45°,  and  ±90°.  The  maximum  error  due  to  the  second  order  terms  will 
occur  at  a  latitude  of  ±45°.  Consequently,  AR^^^  is  the  maximum  error 
which  will  be  caused  by  truncating  second  order  terms : 

AR  -  -  (3/8)  R  f^ 
max  N  /  '  g 
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For  Earth  (R^-6378  km,  f= . 003352 [ 11 ] )  the  maximum  magnitude  of  error  due 
to  second  order  terras  is  27  m.  For  Mars  (R^*=3393  km,  f-«.005185[ll] )  the 
maximum  magnitude  error  due  to  second  order  terms  is  34  m. 

Ellipsoidal  Altitude 

In  this  section,  the  ellipsoidal  altitude  (which  is  the  difference 
between  the  satellite's  radius  and  the  planet's  radius)  will  be 
developed  in  terms  of  the  classical  orbital  elements.  The  radius  of  the 
satellite  is  know  in  terms  of  the  classical  orbital  elements  and  the 
radius  of  the  planet  is  expressed  in  terms  of  the  latitude. 

Equation  (42)  presents  the  equation  for  the  ellipsoidal  altitude 

h  =  r  -  R  (42) 

where  r  is  the  satellite  radius,  and  R  is  the  radius  of  the  planet.  For 
accurate  evaluation  equation  (37)  will  be  employed  to  determine  R  and  r 
will  be  determine  by  numerically  integrating  the  satellite's  equations 
of  motion. 

First  Order  Approximation.  An  approximation  to  the 
ellipsoidal  altitude  may  be  found  by  considering  two-body  motion 
(equation  (3))  and  using  only  first  order  terms  of  equation  (41). 

h  =  a(l  -  e^)/(l  +  e  cosi/)  -  Rg(l  -  f  sin^^)  (43) 

It's  much  more  convenient  to  know  the  ellipsoidal  altitude  as  a  function 

2 

of  orbital  elements  instead  of  latitude;  consequently,  the  sin  ^  term 
will  be  rewritten  in  terms  of  classical  orbital  elements.  Consider  the 
orbit  presented  in  Figure  3  with  inclination  (i),  argument  of  periapsis 
(w) ,  and  true  anomaly  (v).  Using  the  law  of  sines  from  spherical 
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trigonometry,  results  in  an  expression  for  latitude  in  terms  of 
inclination,  argument  of  periapsis,  and  true  anomaly. 

sin<^/sini  =  sin(aH-i/)/sin90 
sinfli  =  sini  sin(c^t/) 


Figure  3.  Relationship  between  Latitude  and  Orbital  Elements 

The  ellipsoidal  altitude  may  now  be  written  in  terms  of  orbital 
elements . 

h  -=  [a(l-e^)/(l  -  e  cosi/)]  -  *  f  sin^i  sin^fuH-i/1]  (44) 

Change  in  Ellivsoidal  Altitude 

The  derivative  of  ellipsoidal  altitude  will  be  zero  for  the 
duration  of  the  constant  altitude  arc.  In  this  section  the  derivative 
of  the  ellipsoidal  altitude  is  calculated  with  respect  to  true  anomaly. 

Differentiating  Equation  (44)  with  respect  to  the  true  anomaly  and 
employing 

sia?a  -  2  cosa  sina 
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results  in  an  expression  for  the  change  in  ellipsoidal  altitude. 


2  2 
dh/di/  =  ae(l-e)  sin»//(l  +  e  cosi/) 

2 

+  f  sin  i  sin[2 (iM-i/) ] 


(45) 


An  arc  of  constant  altitude  would  occur  between  and  when 


It  is  clear  the  altitude  will  remain  constant  throughout  the  entire 
orbit  if  i=e=0.  Also,  as  long  as  the  orbital  elements  remain  constant, 
the  semi -major  axis  only  changes  the  altitude  (not  the  position  or 
length  of  the  arc  of  constant  altitude).  The  argument  of  periapsis  and 
inclination  will  cause  the  position  of  the  constant  arc  to  move  while 
the  eccentricity  must  be  matched  with  the  argument  of  periapsis  and 
inclination  to  produce  an  arc  of  constant  altitude. 

Relationship  between  Mean  and  Osculating  Orbital  Elements 

All  orbital  elements  presented  are  mean  values  with  respect  to 
short  term  variations  caused  by  J2-  When  the  equations  of  motion  are 
numerically  integrated,  the  osculating  elements  are  output  at  each  time 
step  [18].  To  determine  the  secular  change  in  orbital  elements,  the 
mean  orbital  elements  are  needed  instead  of  the  osculating  elements. 

The  method  Kwok  [11],  employs  to  transform  from  mean  to  osculating 
elements  is  presented: 


X.  -  X. '  +  «X. (X. ' ) 
i  1  11 


(46) 
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X^'  are  six  mean  orbital  elements  describing  two-body  motion  and  are 

the  variations  in  orbital  elements  due  to  X^  are  the  osculating 

2 

elements.  Recalling  r  =  p/(l+ecosi/)  .and  p=a(l-e  )  the  change 
in  orbital  elements  may  be  calculated  as  functions  of  the  mean  orbital 
elements . 

6i  =  (3/^)J2(Rg/p' sini'  cosi ' 

X  [  cos2  Cl' '+a; ' )  +  e'  cos(i'  '+2u' )  +  (e’ /3)cos  (3i/ '+2u' )  ] 

6p  =  2  p'  tani '  6i 

60  =  -(3/2)J2(Rg/p' cosi'  [w'  -  M'  +  e'sini/'  -  (l/2)sin2Ci' '-ko') 

-  (e '/2)sinCi^'+3w')  -  (e’ /S)sir\(3i/ '+5ui' )  ] 

6r  =  {(3cos^i'-l) 

X  [2r'(l-e^)/p'  +  e'cosi/ V(l+(l-e'^))  +  1]  -  sin^i'  cos2  (i/ '+w' )} 

6f  =  -(p'/^)  •J2^^e'^^* 

X  [(l-e'2)^/^  +  (pVr')V(l  +  d-e'^)^/^] 

-  2  sin^i'(l  +  e '  cosi/ '  )^sin2  Cv '+w 'd 
Let  F'  =  (l-3cos^I ')/[!  +  (1-e'^)^/^] 

6(i/+w)  =  (J2/8)(Rg/p')^[6(l-5cos^i')(i/'-M') 

+  4e '  siru/ '  (l-6cos^i  '  +  F'  )  +  F'e'^sin2i/' 

+  2(5cos^i  ' -2 )e ' sinCv '+2w'l  +  (7cos^i  ' -l)sin2Ci' '+w '  1 
+  2cos^i 'e '  sinC-3v '+2w '1  ] 
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p  =  p'  +  5p' 


r  =  r'  +  Sr' 
r  =  r'  +  St' 

(l/+{j)  =  (t/'+u')  +  Sil^'+Lo') 

A  -  (p/r)  -  1 

B  =  r 

The  osculating  elements  may  now  be  calculated: 


+  6i' 


+  6Q' 


u  =  tan  ^(B/A) 
a  =  p/(l-e^) 


W  ”  (w+l/)  -  u 


Transformation  from  osculating  to  mean  orbital  elements  requires 
iteration.  First  calculate  the  variations  as  a  function  of  the 
osculating  elements  then  solve  for  the  mean  elements.  Use  these  new 
mean  elements  to  calculate  the  variations.  Continue  iterating  until  a 
group  of  mean  elements  is  found  which  satisfy  equation  (46)  by 
producing  the  original  osculating  elements. 

Impulsive  Velocity  Changes  Needed  to  Correct  Orbit 

Once  an  orbit  with  a  constant  altitude  arc  which  covers  the  correct 
latitude  range  is  chosen  for  a  mission,  secular  changes  to  the  orbital 
elements  will  cause  the  location,  duration,  and  altitude  of  the  arc  to 
change.  Consequently,  the  velocity  impulse  needed  to  negate  the  effects 
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of  the  changes  in  orbital  elements  is  presented.  Figure  4  demonstrates 
the  two  possible  cases  where  the  change  in  semi -major  axis  is  either 
positive  or  negative.  Also  shown  on  the  figure  is  the  change  in 
argument  of  periapsis.  In  practice,  the  optimum  starting  and  ending 
true  anomalies  and  eccentricity  of  the  transfer  ellipse  between  two 
elliptical  orbits  must  be  chosen  to  minimize  the  impulsive  velocity 
change [ 3 : 67 - 71 ] .  In  a  Hohmann  transfer,  the  initial  and  final  orbits 
are  circular  and  the  transfer  ellipse  ranges  from  apoapsis  to  periapsis. 
In  view  of  the  fact  that  the  orbits  being  considered  are  nearly 
circular,  the  portion  of  the  transfer  ellipse  used  in  this  analysis  will 
arbitrarily  be  fixed  to  range  from  apoapsis  to  periapsis.  When  the 
ratio  of  larger  semi-major  axis  to  smaller  semi-major  axis  of  the 
initial  and  final  orbits  being  considered  is  between  0  and  11.94  [3:62], 
the  Hohmann  transfer  is  the  transfer  requiring  a  minimum  velocity 
impulse.  The  cases  to  correct  for  secular  changes  in  orbital  elements 
result  in  a  ratio  of  semi -major  axes  of  about  unity,  with  nearly 
circular  initial  and  final  orbits.  As  a  result,  the  impulsive  velocity 
change  calculated  should  be  near  the  minimum  required  velocity  change. 
The  method  used  to  determine  the  required  velocity  impulse  is  presented. 
From  conservation  of  energy  in  an  elliptical  orbit,  the  magnitude  of  the 
velocity  at  any  point  on  the  ellipse  may  be  determined  by  equation 
(^7): 


V  -  (M(2/r  -  l/a))^/2  (47) 

Also  for  an  elliptical  orbit,  equation  (48)  may  be  employed  to  calculate 
the  angle  (y)  [14:83].  Angle  7  is  measured  from  the  radius  vector  to 
the  velocity  vector. 
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Figure  4.  Geometry  of  Transfer  Orbits 


2  1/2 

siny  =  (l+ecosi/)/(l  +  e  +  2ecosi/)  '  (48) 

2  1/2 

cosy  -  -esiru//(l  +  e  +  Zecosi/) 


Considering  Figure  4,  point  "A"  represents  the  point  where  the  perturbed 
orbit  is  initially  changed  to  follow  the  transfer  orbit,  and  point  "B" 
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is  the  point  where  the  transfer  orbit  is  changed  to  achieve  the  original 
orbit.  The  semi -major  axis  of  the  perturbed  orbit  is  the  semi -major 
axis  of  the  original  orbit  (a)  plus  the  secular  change  in  the  semi -major 
axis  over  one  period  (Aa)  and  similarly  the  eccentricity  of  the  initial 
and  final  orbits  are  e  and  e+Ae  respectively. 

First,  the  velocity  impulse  needed  at  "A"  is  considered.  Equation 
(3)  is  employed  to  calculate  the  radius  of  the  perturbed  orbit  at  "A" 

(r^  p)  where  the  true  anomaly  at  point  "A"  is  Au>  when  Aa<0,  and  180°+Aw 
when  Aa>0.  Knowing  r^  p,  the  velocity  at  "A"  on  the  perturbed  orbit 

(V  )  may  be  evaluated  by  applying  equation  (47). 

A ,  p 

Considering  point  "B" ,  the  radius  of  the  initial  orbit  at  "B" 

(tg  is  calculated  by  employing  equation  (3)  where  the  true  anomaly  is 
180°  when  Aa<0  and  0°  when  Aa>0;  the  velocity  of  the  initial  orbit  at 
point  "B"  is  calculated  by  employing  equation  (47)  with  the  original 
semi -major  axis. 

The  transfer  semi -major  axis  is  found  in  terms  of  the  radius  of  the 
perturbed  orbit  at  point  "A"  and  the  radius  of  the  original  orbit  at 
point  "B" , 


‘t  -  '■'a.p  + 


The  velocity  on  the  transfer  ellipse  may  be  determined  at  points  "A"  and 
"B"  by  applying  equation  (47),  with  a  and  r  (for  V  )  and  r^  .  (for 
V  ).  The  transfer  ellipse  velocity  and  the  initial  ellipse  velocity 
at  point  "B"  are  colinear.  As  a  result,  the  magnitude  of  the  velocity 
impulse  required  at  point  "B"  is  simply  the  difference  between  the 


desired  velocity  (Vg  ^)  and  the  velocity  of  the  satellite  on  the 

transfer  ellipse  (V„  ). 

B ,  t 
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The  velocities  at  point  "A"  are  not  colinear  due  to  the  change  in 

argument  of  periapsis  --  the  angle  between  V  and  V  is  |90°-7  |  ~ 

A  ,  t  A ,  p 

7'.  Consequently, 


'"'a  -  i  (''A.t  •  ''a.p‘=“‘‘1''>^ 


=  [  V  ^  -  2V  V  sin7  +  V  ^ 
A,t  A,t  A,p  A,p 


where  equation  (48)  is  employed  to  calculate  siny  as  a  function  of  the 
true  anomaly  of  the  perturbed  orbit  at  point  A. 

The  velocity  impulse  needed  to  change  the  plane  of  the  orbit  to 
correct  for  the  secular  change  in  inclination  (Ai)  is  now  considered. 

The  velocity  impulse  required  for  this  maneuver  is  directly  proportional 
to  the  magnitude  of  the  velocity  of  the  satellite;  consequently,  the 
maneuver  should  be  accomplished  at  apoapsis  of  the  largest  orbit. 

If  Aa<0 

AV^i  =  2Vg  . IsinfAi/^ll 

and  if  Aa>0 

AV,,.  =  2V.  IsinCAi/f;! 

Ai  A,p'  '  ' 


The  total  change  in  velocity  due  to  secular  changes  in  semi-major  axis, 
argument  of  periapsis,  eccentricity,  and  inclination  is  the  sum  of  the 
velocity  impulses  required  at  points  "A"  and  "B",  and  the  impulse 
required  to  change  the  plane  of  the  orbit. 


AV  -  AV.  +  AV_  +  AV.. 

A  0  Al 


(49) 
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Correcting  secular  changes  due  to  the  change  in  the  longitude  of 
the  ascending  node  was  not  considered.  The  eccentricity  needed  to 
produce  an  orbit  with  a  constant  altitude  arc  is  a  function  of 
inclination,  semi-major  axis,  argument  of  periapsis,  ellipticity, 

and  equatorial  radius  (see  equation  (52)).  The  longitude  of  the  ascending 
node  (fl)  doesn't  affect  the  orbital  elements  needed  to  accomplish  an 
orbit  which  contains  an  arc  of  constant  altitude. 
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Ill  Numerical  Methods 


Search  Method 

Identify  Orbits  with  Arcs  of  Constant  Altitude .  In  order  to  map 
the  combinations  of  mean  orbital  elements  which  produce  orbits  with  arcs 
of  constant  ellipsoidal  altitude,  the  semi-major  axis  (a),  argument  of 
periapsis  (w) ,  inclination  (i),  and  planet  were  chosen  as  independent 
variables.  Several  cases  were  considered  where  the  planet,  semi -major 
axis,  and  perturbing  accelerations  were  fixed.  The  eccentricity  (e) 
which  produced  arcs  of  constant  ellipsoidal  altitude  was  found  as  a 
function  of  argument  of  periapsis  and  inclination.  The  following  matrix 
of  independent  variables  was  considered  for  each  case; 

w  =  j(10)  0<j<18 

i  =  k(10)  0sk<9 


The  method  employed  to  find  values  of  eccentricity  required  an  input  of 
the  ellipsoidal  altitude  as  a  function  of  an  independent  variable  (time 
or  true  anomaly)  over  one  period.  The  change  in  altitude  with  time  (or 
true  anomaly)  was  calculated  at  each  time  step.  The  longest  group  of 
points  where  dh/dt<£  was  found;  where  e  is  the  tolerance  used  to  sort 
the  derivatives.  A  function  (F)  indicating  the  shape  of  this  group  of 
points  was  employed  to  evaluate  the  orbit  and  iterated  to  find  the 
optimum  eccentricity 


m 


m 


n-0 
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where  Ah  is  the  difference  between  the  maximum  and  minimum  altitudes 
max 

attained  in  the  range,  n  is  the  number  of  points  with  dh/dt  -  0,  is 
the  sign  of  the  derivative  at  the  first  point  on  the  range  and  m  is  the 
number  of  times  dh/dt  changes  sign.  The  combination  of  orbital  elements 
which  forced  F  to  be  zero  were  found  (see  Appendix  A  for  a  flow  chart  of 
the  method) . 

Figures  5  and  6  demonstrate  changing  the  eccentricity  changes  the 
satellite's  altitude  above  the  Earth.  For  both  cases  inclination  is 
63.44°  and  semi -major  axis  is  7000  km;  however,  Figure  5  has  an  argument 
of  perigee  of  0°  while  Figure  6  has  an  argument  of  perigee  of  135°. 


i 

I 


Earth:  a  •  7000  kn:  u  =  0  deg;  i  °  63.  44  deg;  ecal  0^  4. 063e~3 
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Figure  5.  Satellite  Altitude  as  a  Function  of  Eccentricity  (w 


0°) 


To  get  an  indication  of  the  sensitivity  of  the  arc  of  constant  altitude 
to  eccentricity,  a  Martian  case  with  a  large  eccentricity  is  also 
presented.  Figure  7  demonstrates  the  changing  altitude  of  a  Martian 
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satellite  as  a  function  of  eccentricity.  The  parameters  which  were  held 

constant  follow:  a  =  3600  km,  w  =  0°,  and  i  =  90°. 

Figures  5,  6,  and  7,  indicate  the  eccentricity  may  differ  from  the 

-3 

optimum  by  as  much  as  .1x10  and  still  produce  a  reasonably  constant 
arc  . 


Earlfi:  »  «  ?006  kn;  u  «  1 35  d*g;  1  «  63.  HH  d*g;  real  c  «  2. 069»“3 
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Figure  6.  Satellite  Altitude  as  a  Function  of  Eccentricity  (w 


135°) 


The  eccentricity  required  to  produce  an  arc  of  constant  altitude  was 

.3 

identified  to  ±  .001  x  10  The  largest  deviation  from  the  fitted  data 

.3 

due  to  the  curve  fit  was  on  the  order  of  ±.01  x  10  (see  Appendix  F) 

.3 

and  most  calculated  curve  fit  points  were  within  ±  .003  x  10  of  the 
data.  The  largest  deviations  occurred  at  larger  eccentricities  where 
larger  deviations  from  the  optimum  eccentricity  may  be  tolerated  and 
still  produce  a  constant  altitude  arc. 
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n«rs;  a  =  3600  kn;  w  •  0  dag:  i  •  90  dag:  acal c  •  7. 668a-3 


Two  methods  were  used  to  produce  the  values  of  ellipsoidal  altitude 
over  one  period.  For  the  two-body  case,  equations  (44)  and  (45)  were 
simply  evaluated  at  each  time  step.  When  perturbative  effects  were 
considered,  ASAP  [11]  was  used  to  numerically  integrate  the  equations  of 
motion  and  provide  values  of  ellipsoidal  altitude  as  a  function  of  time. 
A  central  differences  scheme  was  then  used  to  estimate  the  derivatives 
as  a  function  of  time. 

Curve  Fit.  A  polynomial  fit  was  employed  to  curve  fit  both  the 
two-body  results  and  the  change  in  eccentricity  due  to  J2 .  The  general 
curve  fit  method  is  presented  here  --  independent  variables  which  were 
employed  and  the  coefficients  which  resulted  are  presented  in  Section 
IV. 
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where 


kx  ^ 
k=0  i=0 


C  = 


H: 


k»<0  and  i^O 
5  k=0  or  i=0 

25  k=0  and  i=0 


kx  =  order  of  the  polynomial  in  x 
iy  =  order  of  the  polynomial  in  y 


X  = 


Y  = 


2x  -  (x  +x  .  ) 
_ max  min 

X  •  X 

max  min 


2y  -  (y  +y  •  ) 
-'max  -'min 

^max  ^min 


The  fit  may  only  be  used  for  interpolation,  and  not  extrapolation; 
consequently , 


X  .  <  X  <x 
min  max 


V  •  ^  V  :£y 

mi  J  J  n 


min  •'  max 

The  standard  recurrence  relationship  utilized  to  calculate  the  Chebyshev 
polynomials  is  presented  below 

TqCx)  =  1 
Tj^(x)  -  X 

^n<">  "  2xT^-1<’'>  -  ^n-2^"> 


Choice  of  Independent  Variables .  The  independent  variables  were 
selected  so  a  single  empirical  equation  (good  for  either  planet,  and  all 
the  semi-major  axis  and  Inclination  combinations  considered)  could  be 
developed  to  predict  the  eccentricity  which  produces  orbits  with  arcs  of 
constant  altitude.  The  rational  employed  to  choose  the  non-dimensional 
groupings  of  Independent  variables  for  the  curve  fits  is  now  presented. 
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First,  the  first  order  equation  for  ellipsoidal  altitude  was  considered. 
Equation  (44)  was  rewritten  in  a  non-dimensional  form; 


h  +  R 
_ e_ 

2 

R  fsin  i 
e 


^  ^  1  N-' 


2 

R  fsin  i  '■(l+ecosi/)-' 


2 

+  sin 


h  +  R  =  a 
e 


if'O 


i-0 


It  is  clear  from  this  arrangement  that  for  a  given  inclination  and 

planet,  the  value  of  eccentricity  which  produces  an  arc  of  constant 

2 

altitude  will  depend  on  w  and  (R^fsin  i/a),  while  the  position  of  the 

2 

arc  will  depend  on  i/.  Consequently,  u  and  (R^fsin  i/a)  were  chosen  as 

the  independent  variables  used  to  curve  fit  the  two-body  results  from 

2 

all  the  cases.  In  this  term,  R^f  describes  the  planet,  while  sin  i/a 
describes  the  orbit.  The  change  in  the  optimum  eccentricity  from 
two -body  motion  due  to  the  J2  term  of  the  geopotential  must  also  be 
curve  fit.  The  perturbing  potential  caused  by  the  J2  term  of  the 
geopotential  was  used  to  select  the  portion  of  the  independent  variable 
dependent  on  planet  properties.  V2Q  is  presented  below  for 
completeness . 


20 


Osin^isin^  ((*>+»/)  -  1) 


r 

The  second  zonal  harmonic  and  equatorial  radius  were  chosen  to  represent 
2 

the  planet  (>^2^6  order  to  make  the  variable  non-dimensional, 

2 

the  term  (sini/a)  was  used  to  describe  the  orbit.  The  value  used  as 

2 

the  Independent  variable  was  J2(RgSini/a)  .  Again,  this  term  and  u  were 

used  as  Independent  variables  to  curve  fit  all  the  results  for  the 

change  in  eccentricity  from  the  two-body  cases.  A  correction  for 

was  not  Included.  For  the  Earth  cases  considered,  the  change  in  e 

.3 

caused  by  was  insignificant  (on  the  order  of  .001  x  10  ).  The 


largest  deviation  from  the  J2  results  due  to  occurred  for  Mars 
(a=3600km,  i=90  )  where  the  change  due  to  the  term  was  -.026  x  10 
As  Figure  7  indicates,  a  change  of  this  magnitude  is  not  enough  to 
significantly  change  the  results. 

Cases.  Table  1  presents  the  cases  which  were  considered.  In  each 

case,  the  eccentricity  was  found  which  produced  an  orbit  with  an  arc  of 

constant  altitude  for  the  following  values  of  argument  of  periapsis  and 

inclination:  w  =  j(10)  0<j<18  ,  and  i  =  k(10)  0<k<9 

Without  a  specific  application  to  determine  spacecraft  geometry  and 

attitude,  choosing  a  drag  coefficient,  area,  and  mass  become  arbitrary. 

-  6  2 

For  the  cases  in  Table  1  which  considered  drag,  S/m  =  .02x10  km  /kg, 

-  6  2 

and  Cp  =  2.0  [11]  t/ere  chosen  resulting  in  ^  =  .02x10  km  /kg.  This  is 

a  realistic  estimate  of  as  typical  area  to  mass  ratios  of  satellites 

“6  -  6  2 

(not  including  balloons),  are  between  .002x10  and  .02x10  km  /kg 
[8:xiii] . 

Table  1.  Cases  Considered  to  Determine  e(u),i)  which  Produce 
Constant  Altitude  Arcs 


Case 

Planet 

Semi-Major  Axis  (km) 

Model 

1 

Earth 

6700 

Ivo-Body 

2 

7000 

3 

7500 

U 

6700 

J2 

5 

J3 

6 

J2  +  Drag 

7 

J3  +  Drag 

8 

7500 

J2 

9 

J3 

10 

J2  +  Moon 

11 

J3  +  Moon 

12 

Mars 

3600 

Two -Body 

13 

4000 

14 

3600 

J2 

15 

4000 
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I 

Secular  Considerations 

Calculate  Secular  Changes  in  Mean  Orbital  Elements  over  One  Period. 
Once  the  orbital  elements  which  produce  arcs  of  constant  altitude  were 
identified,  the  secular  changes  in  semi -major  axis,  argument  of 
periapsis,  inclination,  eccentricity  and  longitude  of  the  ascending  node 
were  determined  for  several  cases  by  varying  the  semi-major  axis, 
ballistic  coefficient,  and  planet.  The  start  time  was  not  varied, 
consequently,  cases  for  the  various  positions  of  the  third  body  were  not 
considered.  Perturbations  due  to  zonal  harmonics  up  to  the  sixth  order, 
drag,  and  third  bodies  (Earth:  Moon  and  Sun;  Mars:  Sun)  were  considered. 
For  perturbations  due  to  the  geopotential,  equations  (10)  and  (11)  were 
employed.  The  orbits  are  circular  at  i-=0,  and  equations  (10)  and  (11) 
are  not  valid  for  eccentricities  or  inclinations  of  zero;  consequently, 
the  comparisons  considered  only  the  inclination  range  from  15°  to  90°. 
For  perturbations  due  to  drag,  equations  (30)  through  (34)  were 
numerically  integrated  while  employing  equation  (22)  to  calculate  the 
density.  The  altitude  of  the  satellite  was  calculated  by  applying 
equation  (36),  and  the  rotation  rate  of  the  atmosphere  was  assumed  to  be 
the  same  as  the  planet.  Values  employed  to  calculate  secular  changes  in 
orbital  elements  due  to  drag  and  geopotential  are  presented  in  Table  2. 
For  third  body  effects,  equations  (16)  through  (20)  were  employed.  This 
model  requires  the  position  of  the  third  body  be  provided  with  respect 
to  an  inertial  reference  frame  centered  on  the  planet  with  the  x  axis 
aligned  with  the  vernal  equinox.  The  position  of  the  third  body  was 
calculated  once  per  orbit.  The  position  was  found  as  a  function  of  time 
(in  Julian  centuries  since  epoch  (T)).  The  time  used  to  define  the 
position  of  the  third  body  was  a  start  time  (input  variable)  plus  half 
the  period  calculated  using  the  original  mean  orbital  elements. 
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Table  2.  Values  used  to  Calculate  Secular  Changes  in  Orbital  Elements 


Parameter 

Earth 

Mars 

Reference 

j2 

? 

? 

,5 

"6 

1.082627  X  10‘® 

-2.542  X  10'^ 

-1.609  X  10'^ 

-2.192  X  10‘5 

5.23  X  lO'^ 

1.960454  X  10'® 
-3.144926  X  10'^ 
-1.889437  X  10'^ 
2.669249  x  lO'^ 
-1.340757  X  10*^ 

Earth:  [20] 
GEM- 6  Model 

Mars :  [ 3 ] 

H  (km)  - 

Pn  ) 

Zq  (km) 

78.4 

5.215  X  10 

500 

9.98 

1.72  X  10 

100 

Earth:  [11] 
Mars  :  [10] 

u>  (rad/s) 

7.292  X  10'^ 

7.088  X  10"^ 

[11] 

The  positions  which  were  used  are  presented. 

Position  of  the  Sun  with  respect  to  the  Earth  (ref  JD:  2415020 .)[ 5 ] : 

a  =  1.496  X  10®  km 
s 

e  =  1.675104e-2  -  1.1444e-5T  -  9.4e-9T^ 
s 

i  =  .4093197474  -  6.217910e-5T  -  2.1468e-9T^  +  1.7977e-10T®  rad 
s 

0=0  rad 
s 

w  =  4.908229653  +  8 . 24149855e- 7T  +  5.9167e-7T^  +  1.22e-9T®  rad 
s 

M  =  6.25658  +  172.0197T  -  1.9548T^  -  1.22e-9T^  rad 
s 


Position  of  the  Moon  with  respect  to  the  ecliptic  --  all  angles  are  in 
radians,  (ref  JD;  2415020.) 

a  =  384.4  X  10  ®  km 
m 

e  -  .054900489 
m 

i  '  -  8.980411e-2 
m 

0  '  -  4.523601514  -  9.242202945T  +  2 . 7174776e-6T^  +  8.7266e-10T® 
m 

w  '  -  5.8351515  +  19.44368005T  -  1.35071e-5T^  -  4.53786e-9T^  -  0  ' 
m  m 

M  '  -  4.71997  +  2299. 715T  -  1.4835e-6T^  +  6.81e-10T®  -  Q  '  -  w  ' 
m  mm 
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The  position  of  the  Moon  with  respect  to  the  Earth  is  required. 

Transformations  to  determine  the  position  of  the  Moon  with  respect  to 

the  Earth  instead  of  the  ecliptic  are  completed  using  spherical 

trigonometry  and  are  presented  below.  Figure  8  demonstrates  the 

geometry  of  the  orbits [ 5 : 107- 109 ] .  To  find  i^  use  spherical 

trigonometry  to  find  a, 

cosQ  =  -cosi  cosi  +  sini  sini  '  cosfl 
s  m  s  m  m 

where  a  is  in  the  second  quadrant.  Now  i^  may  be  found, 


i  =  tr  —  a 
m 


Ecliptic  Plane 


n  n 


Plan*  of  Earih's  Equaior 


U*mal  Equi  nox 


Moon'*  Orbi  i 


Figure  8.  Geometry  of  Lunar  Orbit  with  Respect  to  Ecliptic  and  Earth 
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In  order  to  find  solve  for  "b"  using  the  law  of  cosines  for  angles 

and  ensure  "b”  lies  in  the  same  hemisphere  as  C2  ' 

^  m 

cash  =  (cosi  +  cosa  cosi  ')/(sina  sini  ') 
s  w  m 

w  -  w  '  +  b 
m  m 

Finally,  solve  for  using  the  law  of  cosines  for  sides  and  ensure 
n  is  in  the  same  hemisphere  as  0  ' . 

cosCI  =  cosQ  '  cost  +  sinQ  '  sinb  cosi  ' 
mm  m  m 

Position  of  the  Sun  with  respect  to  Mars  (reference  Julian  date 

2448400. 5) [11] 

a  =  227.941  x  10^  km 
s 

e  =  .933969  x  10*^ 
s 

ig  -  25.191153  deg 

=  0.  deg 

=  -109.05076  deg 
=  -171.60476  deg 

The  equations  defining  perturbations  due  to  the  third  body  require 
the  semi-major  axis,  radius,  right  ascension,  declination  and  mean 
motion  of  the  third  body.  The  methods  which  were  employed  to  find  these 
values  are  presented.  The  radius  is  calculated  using  equation  (3)  where 
1/  is  approximated  by  an  expansion  in  mean  anomaly  and  eccentricity [ 15 ] 

1/  =  M  +  (2e  -  eV^)sinM  +  (5eV^  +  lle^/2U)sin2H 

+  (13/12)e^sin3/f  +  (103/96)e^sin4/f  (50) 
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The  mean  motion  is  found  using  the  semi -major  axis  and  mass  of  the  third 
body: 


/  /  3,1/2 
n  =  (/i/a  ) 


The  declination  ($')  and  right  ascension  (O')  are  calculated  using 
relationships  from  spherical  trigonometry.  The  range  of  inclinations  of 
the  third  bodies  considered  in  this  thesis  were  only  in  the  first 
quadrant;  consequently,  declination  is  in  the  same  hemisphere  as  (uH-i/) . 
The  right  ascension  is  the  sum  of  the  longitude  of  the  ascending  node 
(0)  and  "b" .  Equation  (51)  is  employed  to  calculate  "b"  which  must  be 
in  the  same  quadrant  as  (uf+u) . 

=  sin  ^[sini  sin(w+»/)] 
tanb  =  cosi  tan$'  (51) 


Determine  Orbit  Stability.  Finally,  in  order  to  compare  the 
various  orbits  to  determine  which  was  the  most  stable,  the  change  in 
velocity  needed  to  transfer  the  satellite  from  the  perturbed  orbit  back 
to  the  original  orbit  was  considered.  Equation  (49)  was  used  to 
determine  the  change  in  velocity. 
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Cases.  Table  3  presents  the  cases  consider  to  determine  secular 


changes  in  orbital  elements. 


Table  3.  Cases  Considered  to  Determine  Secular  Changes  in  Orbital 
Elements  and  Velocity  Change  as  a  Function  of  w  and  i 
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IV  Results 


Two -Body  Solutions 

First,  the  eccentricity  which  produces  constant  altitude  arcs  was 
found  as  a  function  of  argument  of  periapsis  and  inclination  while 
holding  the  semi-major  axis  and  planet  fixed.  Figures  B-1  and  B-2  present 
the  results  for  Earth  with  semi -major  axes  of  6700  km  and  7500  km 
respectively.  Figures  B-3  and  B-4  present  the  results  for  Mars  with 
semi -major  axes  of  3600  km  and  4000  km. 

The  eccentricities  identified  are  all  quite  small  and  the  orbits 
are  almost  circular.  As  may  be  seen  in  all  the  figures,  Se/Sus  is 
greatest  at  argument  of  periapses  of  0°  and  180*^  and  a  minimum  at  45°, 
and  135°;  and  6i/6w  is  smallest  for  inclinations  of  0°  and  90°.  As 
would  be  expected  from  equation  (44),  the  eccentricity  increases  as  the 
inclination  increases,  and  is  zero  when  the  inclination  is  zero. 

Comparing  Figures  B-1  and  B-2,  and  Figures  B-3  and  B-4,  indicates  the 
eccentricity  decreases  as  the  semi-major  axis  increases. 

Perturbed  SoluCion 

Figures  B-5  though  B-8  present  the  change  in  eccentricity  due  to 

Although  all  the  perturbations  listed  in  Table  1  were  considered, 

it  became  apparent  for  a  single  period,  the  only  significant  effect  on 

the  eccentricity  required  to  achieve  an  arc  of  constant  altitude  was  due 

to  the  J2  term  of  the  geopotential.  For  Earth,  the  difference  between 

the  eccentricity  found  considering  J^,  drag,  and  the  Moon  and  the 

.3 

eccentricity  found  considering  J2  was  on  the  order  of  .001  x  10  For 

Mars,  the  largest  change  in  eccentricity  due  to  was  on  the  order  of 

.3 

.01  X  10  .  Changes  due  to  drag  weren't  considered  for  Martian  cases. 

The  Martian  atmosphere  is  much  less  dense  than  the  Earth's  and 
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short  term  perturbations  due  to  drag  shouldn’t  significantly  affect  the 
choice  of  eccentricity. 

The  change  in  eccentricity  presented  in  Figures  B-5  through  B-8  is 
subtracted  from  the  eccentricity  predicted  using  the  two-body  motion  to 
give  the  required  eccentricity.  Figures  B-5  and  B-6  present  the  results  for 
Earth  with  semi-major  axes  of  6700  km  and  7500  km  respectively.  Figures 
B-7  and  B-8  present  the  results  for  Mars  with  semi -major  axes  of  3600  km 
and  4000  km . 

Once  the  optimum  eccentricity  was  found  which  produced  constant 
altitude  arcs,  the  location,  duration,  and  latitude  range  covered  by  the 
arcs  were  found.  The  mid-latitude  presented  is  the  average  of  the 
maximum  and  minimum  latitudes  covered  by  the  constant  altitude  arcs,  while 
the  latitude  range  presented  is  the  difference  between  the  maximum  and 
minimum  latitudes  covered.  Both  the  duration  of  the  arc  and  the 
latitude  range  of  the  arc  are  presented  for  a  maximum  change  in  altitude 
of  100  m  and  1  km.  The  change  in  altitude  is  the  difference  between  the 
maximum  and  minimum  altitudes  attained  by  the  arc.  Figures  C-1  through 
C-12  present  the  results  for  Earth  where  Figures  C-1  through  C-6  are  for 
the  6700  km  semi-major  axis  orbits  and  Figures  C-7  through  C-12  are  for 
the  7500  km  orbits.  Figures  C-13  through  C-24  are  the  results  for  the 
Martian  cases  where  Figures  C-13  through  C-18  present  the  3600  km 
results  and  Figures  C-19  through  C-24  are  for  the  4000  km  semi-major 
axis  case. 

Mid -Latitude  of  Arcs.  Figures  C-1,  C-7,  C-13,  and  C-19  present  the 
mid-latitude  as  a  function  of  argument  of  periapsis  and  inclination. 

There  are  a  large  number  of  combinations  which  produce  arcs  with  a  given 
mid-latitude.  The  arcs  are  positioned  over  the  northern  hemisphere  when 
the  argument  of  periapsis  is  between  0°  and  180°  and  are  positioned  over 
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the  southern  hemisphere  when  the  argument  of  periapsis  is  between  180° 
and  360°.  As  inclination  increases,  the  magnitude  of  the  mid-latitude 
of  the  arc  increases.  Of  course,  the  inclination  must  be  at  least  as 
large  as  the  sum  of  the  mid-latitude  and  half  the  desired  latitude 
range.  The  mid-latitude  was  calculated  by  taking  the  difference  between 
the  maximum  and  minimum  latitudes  achieved  by  the  arc;  consequently,  the 
mid-latitude  of  polar  orbits  with  the  arc  positioned  over  the  poles  (w  = 
90°  and  270°)  is  less  than  90°. 

Altitude  of  Arc.  Figures  C-2,  C-8,  C-14,  and  C-20  present  the 
altitude  of  the  arcs  as  a  function  of  argument  of  periapsis  and 
inclination.  The  altitude  is  dependent  on  the  mid-latitude  of  the  arc, 
the  argument  of  periapsis ,  and  the  true  anomaly  range  covered  by  the 
arc.  At  argument  of  periapses  of  0°  and  180°,  the  constant  altitude  arc 
occurs  at  apoapsis,  while  the  arc  occurs  at  periapsis  when  the  argument 
of  periapses  are  90°  and  270°. 

Duration  of  Arc.  Figures  C-3,  C-9,  C-15,  and  C-21  present  the 
duration  of  the  arcs  where  the  maximum  change  in  altitude  is  less  than 
100  m.  When  i“e=0,  the  arc  lasts  the  entire  orbital  period;  however, 
the  latitude  range  covered  is  zero.  The  duration  of  the  arc  is 
reasonably  constant  for  inclinations  larger  than  20°  and  for  a  given 
argument  of  periapsis  the  duration  decreases  as  the  inclination 
increases.  Longer  duration  arcs  may  be  produced  when  the  argument  of 
periapsis  is  an  integer  multiple  of  90°  with  the  largest  duration  arcs 
occurring  at  w  -  0°  or  180°.  When  the  argument  of  periapsis  is  0°  or 
180°,  the  mid-latitude  of  the  constant  altitude  arc  is  positioned  at  0° 
latitude.  For  these  orbits,  the  symmetry  of  the  orbit  and  planet  are 
aligned  causing  a  longer  constant  altitude  arc.  Figures  C-4,  C-10, 

C-16,  and  C-22  present  the  duration  of  the  arcs  when  the  change  in 
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altitude  is  less  than  1  km.  As  the  allowable  altitude  variation  of  the 
arc  is  increased,  the  duration  increases. 

Latitude  Range  of  the  Arc.  Figures  C-5,  C-11,  C-17,  and  C-23 
present  the  latitude  range  of  the  arcs  as  a  function  of  argument  of 
periapsis  and  inclination  when  the  altitude  variation  is  less  than  100 
m.  Figures  C-6,  C-12,  C-18,  and  C-24  present  the  latitude  range  when 
the  altitude  variation  is  less  than  1  km.  In  all  cases,  the  largest 
latitude  coverage  occurs  when  the  argument  of  periapsis  is  0°  or  180°. 
Again,  this  is  due  to  the  matched  symmetry  of  the  orbit  and  planet,  and 
coincides  with  the  longest  duration  arcs.  Also,  at  a  given  argument  of 
periapsis,  the  latitude  coverage  increases  with  increasing  inclination. 
The  latitude  range  presented  for  an  inclination  of  90°  at  argument  of 
periapses  of  90°  and  270°  is  low  due  to  the  method  used  to  calculated 
the  latitude  range.  An  orbit  with  these  parameters  will  have  the 
constant  altitude  arc  positioned  over  the  poles  and  when  the  difference 
between  the  maximum  and  minimum  latitudes  achieved  is  calculated,  the 
resulting  latitude  range  will  be  less  than  the  range  reported  if  the  arc 
was  not  positioned  over  a  pole.  A  latitude  range  of  about  50°  is 
possible  while  maintaining  an  altitude  variation  of  less  than  100  m  and 
a  latitude  range  of  up  to  90°  can  be  achieved  if  an  altitude  variation 
of  less  than  1  km  is  acceptable. 

Non-dimensional ized  Results 

The  results  from  the  previous  cases  were  combined.  Figure  9 

2 

presents  a  plot  of  eccentricity  as  a  function  of  u>  and  R^fsin  i/a  for 

the  two-body  cases,  while  Figure  10  presents  the  change  in  eccentricity 

2 

from  the  two-body  cases  as  a  function  of  u,  and  J2(RgSini/a)  . 

The  results  presented  in  Figures  9  and  10  were  curve  fit  producing 
a  single  empirical  equation.  Equation  (52)  may  be  used  with  the 
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coefficients  presented  in  Tables  4  and  5  to  determine  the  eccentricity 
required  to  produce  an  arc  of  constant  altitude  given  the  semi -major 
axis,  argument  of  periapsis,  inclination,  ellipticity,  and  the  J2  zonal 
harmonic  of  the  geopotential . 


e  =  -  Ae  i#0  (52) 

e  =  0  i=0 


The  inclination  is  combined  with  planet  properties  and  differs  for  e2g 
and  Ae.  For  ©28’  ^2B  ^  function  of  semi -major  axis,  inclination, 

ellipticity,  and  equatorial  radius 
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1.0 


2B 


^^2B^^2B(max) 


2  . 


2B 


(f  Rg/a)  sin  i 


0  <  i.„  <  .48874383  x  10 
Zd 


-2 


2B(max) 


For  Ae ,  1^2  is  a  function  of  semi-major  axis,  inclination,  ellipticity, 
and  J2 


0  < 


I 


J2 


J2  ^^J2/^J2(max)  ' 

‘j2  ■  •’2  <V®>^ 

<  .17418950  X  10-2  . 


Table  4.  Coefficients  for  e 


\,i 

i 

0 

1 

2 

k 

0 

1.3382318E-2 

6.6831547E-3 

-8.8992362E-6 

1 

5.9572268E-5 

3.9537390E-5 

9.6394148E-6 

2 

2.2536540E-3 

1.1263017E-3 

-2.0021737E-6 

3 

4.6696319E-5 

3.0815889E-5 

7.3419897E-6 

4 

5.9702036E-4 

3.0177004E-4 

3.1364165E-6 

5 

2.6463970E-5 

1.7652256E-5 

4.3897306E-6 

6 

1.8920545E-3 

9.6492837E-5 

2.2358554E-6 

7 

8.4988295E-6 

5.9054707E-6 

1.6974386E-6 

8 

4.2832099E-5 

2.2063680E-5 

1.2341417E-6 
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Table  5.  Coefficients  for  Ae 


2.4192677E-3 

2.2288511E-5 

4.6494823E-4 

1.7652175E-5 

1.7222103E-4 

9.8102319E-6 

6.7878651E-5 

2.9070219E-6 

1.8819281E-5 


1.2121056E-3 

1.4267904E-5 

2.2535319E-4 

1.1116183E-5 

8.6194009E-5 

6.3723936E-6 

3.4179548E-5 

2.1377139E-6 

8.8447405E-6 


1.0566496E-5 

2.9502728E-6 

5.4081346E-6 

2.2108308E-6 

6.6663719E-6 

1.3617805E-6 

3.6139002E-6 

5.7416986E-7 

9.3013404E-7 


1.3977998E-5 

2.0193706E-8 

1.2653309E-5 

1.9007614E-7 

7.1543029E-6 

■7.8643838E-8 

3.1601894E-6 

■2.5954531E-7 

1.0599178E-6 


4.5367581E-8 

1.5042158E-7 

-3.2310241E-7 

8.1730486E-8 

1.7824329E-6 

8.5013670E-8 

1.0169085E-6 

7.5947404E-8 

•2.6786479E-7 


-2.4362035E-6 

-1.5487548E-7 

9.7169557E-7 

-4.1539688E-7 

2.3768136E-6 

7.9025722E-8 

2.2980839E-6 

4.2068010E-7 

1.4068201E-6 


Figures  11  and  12  are  plots  of  equations  (53)  and  (54).  Comparing 
Figures  11  and  12  with  Figures  9  and  10  Indicates  the  curve  fit  is 
sufficiently  accurate.  Tabulated  comparisons  of  the  data  points  and 
eccentricities  calculated  using  equation  (52)  are  presented  in  Appendix 
F.  A  discussion  of  the  required  accuracy  in  predicting  eccentricity  is 
presented  in  Section  III. 
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Earth  &  tiara:  aCw.'f*R*(  am  i)"2/a3 
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ARGUTEHT  OF  PERIflPSlS  (dag) 


Figure  9.  Non-dimensional  Presentation  of  Two-Body  Solutions 


Earth  &  tiara:  dalta  aCy.  J2*(  R*(  ai  n  i)/a)'‘Z] 


0.  50.  100.  150.  200. 

0R(»CNT  OF  PERIRPSIS  (dag) 


Figure  10.  Non-dimensional  Presentation  of  Change  Due  to  J., 


OELTfl  ECCEHTRICITV 


Curv*  Fit  (Kvplvr):  •(ii.('f*R/a)*(*iru''2)) 


0.000  - -  - - 0. 0000 

0.  50.  100.  150.  20a 

flRGUHEHT  OF  PERIflPSIS  (d*g) 


Figure  11.  Curve  Fit  of  Two-Body  Solutions 


Corv*  hi  (J2):  •(«.J2*(R*sini/«)''2) 


Figure  12.  Curve  Fit  of  Change  Due  to  J2 


53 


Secular  Considerations 

The  secular  changes  in  orbital  elements  were  calculated  using  the 
mean  orbital  elements  which  produce  orbits  with  arcs  of  constant 
altitude.  Perturbations  due  to  the  geopotential,  drag,  and  third  bodies 
were  considered.  The  secular  changes  in  a,  e,  i,  w,  and  0  are  presented 
as  a  function  of  the  initial  argument  of  periapsis  and  inclination.  In 
all  cases,  the  initial  longitude  of  ascending  node  was  set  to  0°. 

Secular  Change  in  Semi-Major  Axis  over  One  Period.  Figures  D-1 
through  D-6  present  the  results  for  the  secular  change  in  semi-major 
axis.  Figures  D-1  and  D-2  present  the  results  for  Earth  with  a 

-  6  2 

semi-major  axis  of  6700  km  and  ballistic  coefficients  of  .01x10  km  /kg 
•  6  2 

and  .04x10  km  /kg  respectively  while  Figures  D-3  and  D-4  present  the 

same  results  for  Mars  with  a  semi-major  axis  of  3600  km.  In  order  to 

consider  cases  in  which  drag  doesn't  play  an  important  role,  Figures  D-5 

and  D-6  present  the  results  for  Earth  (a  “  7500  km)  and  Mars  (a  =  4500 

“6  2 

km)  with  a  ballistic  coefficient  of  .04x10  km  /kg.  Comparing  Figures 
D-1  and  D-2  and  Figures  D  3  and  D-4  indicate  increasing  the  ballistic 
coefficient  increases  the  change  in  semi-major  axis  for  low  orbits. 
Equation  (30)  is  directly  proportional  to  the  ballistic  coefficient  and 
the  square  of  the  semi -major  axis,  and  supports  this  observation.  For 
both  Earth  and  Mars,  the  minimum  change  in  semi -major  axis  due  to 
atmospheric  drag  occurs  at  an  argument  of  periapsis  of  45°  or  225°  and 
an  inclination  of  63°.  In  Figures  D-5  and  D-6,  the  geopotential  caused 
the  most  significant  change  in  semi-major  axis.  Equation  (10)  indicates 
the  change  in  semi -major  axis  due  to  the  geopotential  will  be  zero  when 
u  is  0°  or  180°  or  i  is  63.44°.  These  areas  where  semi -major  axis 
doesn't  change  are  clearly  presented  by  Figures  D-5  and  D-6.  The  small 
negative  offset  of  Figure  D-5  is  due  to  atmospheric  drag. 
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Secular  Change  in  Inclination  over  One  Period.  Figures  D-7  through 

D-12  present  the  results  for  the  change  in  inclination.  Again,  Figures 

D-7  and  D-8  are  for  Earth  (a=6700kni)  and  Figures  D-9  and  D-10  are  for 

Mars  (a=3600kn)) ,  where  the  first  figure  for  each  planet  presents  the 
•  6  2 

data  for  ^=.01x10  km  /kg  and  the  second  figure  presents  data  for 
-  6  2 

/8=. 04x10  km  /kg.  Figures  D-11  and  D-12  present  results  for  Earth 
(a=7500km)  and  Mars  (a=4500km)  with  /3=. 04x10  ^  km^/kg.  For  Earth  the 
largest  change  in  inclination  was  caused  by  the  Moon  in  Figures  D-7  and 
D-11  and  drag  in  Figure  D-8.  For  Mars  the  largest  change  in  inclination 
was  caused  by  the  geopotential  and  Sun  in  Figure  D-9,  drag  in  Figure 
D-10,  and  geopotential  at  low  inclinations  and  the  Sun  at  high 
inclinations  in  Figure  D-12.  Comparing  Figures  D-7  and  D-8  and  Figures 
D-9  and  D-10,  indicates  at  low  altitude,  high  inclination  and  argument 
of  periapsis,  drag  significantly  affects  the  change  in  inclination  by 
reducing  the  inclination.  The  secular  change  in  inclination  due  to  drag 
is  calculated  using  equation  (32)  which  is  zero  when  the  inclination  is 
zero.  If  all  other  parameters  remain  constant  equation  (32)  will  have 
the  largest  magnitude  when  the  inclination  is  90°. 

Secular  Change  in  Argument  of  Periapsis  over  One  Period.  Figures 
D-13  through  D-18  present  the  results  for  the  secular  change  in  argument 
of  periapsis.  Figures  D-13  and  D-14  are  for  Earth  (a-6700km)  and 
Figures  D-15  and  D-16  are  for  Mars  (a-3600km)  --  for  each  planet,  •= 
.01x10  km  /kg  is  presented  first  while  P  -=  .04x10  km^/kg  is 
presented  second.  Figures  D-17  and  D-18  present  results  for  Earth 
(a-7500km)  and  Mars  (a-4500km)  with  P  -  .04x10  ^  km^/kg.  In  all  cases 
considered,  the  largest  secular  change  in  the  argument  of  periapsis  was 
caused  by  the  geopotential.  By  comparing  the  figures,  it  may  be 
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concluded  that  the  change  in  argument  of  periapsis  is  not  significantly 
affected  by  drag.  The  large  magnitudes  at  low  inclinations  were  from 
the  ,  and  ^22  components  of  equation  (11).  The  change  in  argument 

of  periapsis  is  zero  at  an  inclination  of  63.4°  when  the  argument  of 
periapsis  is  zero.  Mars  and  Earth  have  the  opposite  sign  on  the  change 
in  argument  of  periapsis  at  high  inclinations  --  this  is  due  to  the 
difference  in  signs  of  the  terms  of  their  respective  geopotentials. 

Secular  Change  in  Eccentricity  over  One  Period.  Figures  D-19 
through  D-24  present  the  results  for  the  change  in  eccentricity.  The 
cases  presented  are  in  the  same  order  as  previously  mentioned.  In  the 
cases  with  the  smallest  semi-major  axes  and  largest  ballistic 
coefficient  for  both  Earth  and  Mars  (Figures  D-20  and  D-22),  the  change 
in  eccentricity  was  due  primarily  to  drag.  For  all  the  other  cases 
considered,  the  change  in  eccentricity  was  primarily  due  to  the 
geopotential.  Drag  increases  the  reduction  is  eccentricity  --  Figure 
D-20  has  been  shifted  down  relative  to  Figure  D-19,  and  similarly  Figure 

D-22  has  been  shifted  down  relative  to  Figure  D-21. 

Secular  Change  in  Longitude  of  the  Ascending  Node  over  One  Period. 
Figures  D-25  through  D-30  present  the  changes  in  longitude  of  the 
ascending  node.  The  change  in  longitude  of  the  ascending  node  was 
primarily  due  to  the  geopotential  in  all  the  cases  considered.  All  the 
figures  indicate,  the  minimum  change  occurs  at  an  inclination  of  90°  and 
the  largest  magnitude  change  occurs  at  an  inclination  of  0°.  The  Q2 

term  of  equation  (11)  is  a  multiple  of  cosi,  and  consequently  indicates 

the  secular  change  in  longitude  of  the  ascending  node  will  be  largest 
when  the  inclination  is  0°  and  zero  when  the  inclination  is  90°. 
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Velocity  Impulse  Required  to  Counteract  Secular  Changes.  Figures 
E-1  through  E-6  present  the  velocity  impulses  needed  to  move  the 
satellite  from  the  perturbed  orbit  to  the  original  orbit.  For  all  the 
cases  considered,  the  change  in  velocity  was  driven  by  the  required 
plane  change  to  correct  inclination.  As  a  result,  station  keeping  fuel 
can  be  minimized  by  choosing  an  orbit  which  minimizes  the  change  in 
inclination. 
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V  Conclusions 


In  this  thesis  an  equation  was  developed  to  empirically  determine 
the  mean  orbital  elements  which  produce  arcs  of  minimum  altitude 
variation  about  an  oblate  planet.  These  orbits  are  useful  for 
surveillance  or  scientific  missions  using  optics  with  fixed  focal 
lengths.  The  resulting  orbits  were  considered  and  the  characteristics 
of  minimum  variation  arcs  were  identified.  The  secular  changes  to  the 
mean  orbital  elements  due  to  geopotential,  drag,  and  third  body 
perturbations  were  found  and  the  velocity  impulse  needed  to  maintain  the 
desired  orbit  was  also  calculated. 

Conclusions 

1.  Orbits  which  contain  arcs  of  minimum  altitude  variation  were 
identified . 

2.  An  empirical  relationship  was  developed  to  identify 
combinations  of  the  mean  orbital  elements  which  produce  arcs  of  constant 
altitude.  The  equation  determines  eccentricity  as  a  function  of 
semi-major  axis,  inclination,  argument  of  periapsis,  the  J2  zonal 
harmonic,  and  the  planet's  ellipticity  and  equatorial  radius. 


e  =  -  Ae  i^'O  (52) 

e  =  0  i=0 


where 


'2B 


8  2 

I  I 

k-0  i-0 


^  ^k,i 


(53) 
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8  5 

k=0  i-0 

The  coefficients  for  the  equation  are  presented  in  Tables  4  and  5,  and 
the  independent  variables  are  defined  in  Section  IV,  This  curve  fit  is 
valid  for  inclinations  from  0°  to  90°  for  Earth  orbits  with  semi -major 
axes  greater  than  one  Earth  radius  and  Martian  orbits  with  semi -major 
axes  greater  than  3600  km.  For  smaller  Martian  semi-major  axes,  the 
inclination  range  becomes  limited.  For  Martian  orbits  with  semi-major 
axes  of  one  planet  radius,  the  equation  is  valid  for  an  inclination  range 
of  0°  to  70°.  The  valid  inclination  range  increases  as  the  semi -major 
axis  increases  until  the  full  range  of  0°  to  90°  is  achieved  for 
semi -major  axes  greater  than  or  equal  to  3600  km. 

3.  All  the  orbits  with  arcs  of  minimum  altitude  variation  were 
nearly  circular.  The  orbits  became  more  eccentric  for  smaller 

semi -major  axes  and  larger  inclinations. 

4.  The  velocity  impulse  required  to  counteract  secular  changes 
to  the  mean  orbital  elements  due  to  a  non-spherical  gravitational  field, 
drag,  and  third  bodies  was  calculated.  In  all  cases,  the  magnitude  of 
the  required  impulse  was  driven  by  the  plane  change  to  correct  the 
inclination.  Consequently,  in  order  to  minimize  station  keeping  fuel, 
an  orbit  with  the  minimum  change  in  inclination  should  be  chosen. 

Recommendations 

1 .  Further  work  should  be  done  to  identify  orbits  with  arcs  of 
minimum  altitude  variation  over  specific  latitude  and  longitude  ranges. 

2.  Third  body  perturbations  of  orbits  with  large  semi -major  axes 
should  be  considered  for  different  third  body  positions. 
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SUBROUTINE  FCALC 


Figures  B-1  through  B-8  present  the  optimujD  eccentricity  to  produce 
arcs  of  minimuin  altitude  variation. 

Figures  B-1  and  B-5  present  results  for  Earth  orbits  with  a 
semi-major  axis  of  6700  km;  Figures  B-2  and  B-6  present  results  for 
Earth  orbits  with  7500  km  semi-major  axes;  Figures  B-3  and  B-7  present 
Martian  orbits  with  3600  km  semi-major  axes;  and  Figures  B-4  and  B-8  are 
Martian  results  with  semi-major  axes  of  4000  km. 

Figures  B-1  through  B-4  present  solutions  to  the  two-body  problem, 
while  Figures  B-5  through  B-8  present  the  change  in  optimum  eccentricity 
caused  by  perturbations  due  to  J2 .  The  optimum  eccentricity  for  a  given 
argument  of  periapsis  and  inclination  is  found  by  subtracting  the  change 
due  to  J2  from  the  two-body  solution. 
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Earth  (a  «  7500  kn):  T«o-Body  Solution 


RRGUnEMT  OF  PERinPSIS  (dag) 


I  Figure  B-2.  Earth  (a-7500  km);  Eccentricity  to  Produce  a  Constant  Arc 

I 


64 


rUrs  (a  •  3686  kn):  Two~Bo<ig  Solu'tion 


Figure  B 


Figure 


0.  58.  100.  158.  208. 

ARGUMEHT  OF  PERI  APSIS  (dag) 


-3.  Mars  (a-3600  km);  Eccentricity  to  Produce  a  Constant  Arc 


rUrs  (a  •  H868  kn):  Two-Bodg  Solirhon 


ARGUnEMT  OF  PERIRPSIS  (dag) 


B-4.  Mars  (a-4000  km):  Eccentricity  to  Produce  a  Constant  Arc 
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DELTA  ECCENTRICITV  ,  DELTA  ECCEHTRICITV 


Ear-lh  (a  >  E?00  kn):  DiTfaranca  Ba'Luaan  Two-6odg  and  J2  Soluiions 
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B-5.  Earth  (a=6700  km);  Change  in  Eccentricity  due  to  J2 
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Figure  B-6.  Earth  (a-7500  km):  Change  in  Eccentricity  due  to  J 
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B-7.  Mars  (a=3600  kre) :  Change  in  Eccentricity  due  to  J2 
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B-8.  Mars  (a-4000  km):  Change  in  Eccentricity  due  to  J 
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Appendix  C 

Characteristics  of  Arcs  of  Minimum  Altitude  Variation 

Figures  C-1  through  C-24  present  the  characteristics  of  the  arcs  of 
minimum  altitude  variation. 

Figures  C-1  through  C-6  present  the  results  for  Earth  with  a 
semi-major  axis  of  6700  km.  Figures  C-7  through  C-12  present  results 
for  Earth  with  a  semi -major  axis  of  7500  km.  Figures  C-13  through  C-18 
present  results  for  Mars  with  a  semi-major  axis  of  3600  km  while  Figures 
C-19  through  C-24  are  for  Martian  orbits  with  4000  km  semi-major  axes. 

For  each  of  the  cases  considered,  the  first  figure  presents  the 
mid-latitude  of  the  constant  altitude  arc,  and  the  second  figure 
presents  the  altitude  of  the  arc.  The  third  and  fourth  figures  present 
the  duration  of  the  arc  for  a  change  in  altitude  less  than  100  m  and 
1  km  respectively.  Finally,  the  last  two  figures  present  the  latitude 
range  of  the  arc  for  a  change  in  altitude  of  less  than  100  m  and  1  km. 


68 


niD  LftTITUDE  (deg) 


Eerth  ( •“6?00t<n) :  tti  d  Latitude  ot  ftrc  of  Constant  Altitide 


66. 


-50.  ^  -  d0. 

-  20.  IHCLINftTION  (deg) 

-100.  - . - - - - -  0- 

0.  100.  200.  380.  H00. 

ftRGUnENT  OF  PERI0PSIS  (deg) 

Figure  C-1.  Earth  (a=6700  km):  Mid-Latitude  of  Arc 
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Figure  C-2.  Earth  (a-6700  km):  Altitude  of  Arc 


69 


Ear-Wi  ( asEPBOkn) :  DuraEi  on  of  Consfani  R1  f i  ii  da  ( dal  ta  h  <  .  1  kn) 
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Figure  C-3.  Earth  (a-=6700  km):  Duration  of  Arc  Ah  <  .  1  km 
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Figure  C-4.  Earth  (a-6700  km):  Duration  of  Arc  Ah  <  1km 
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Figure  C-fa.  Earth  (a-6700  km):  Latitude  Range  of  Arc  Ah  <  1km 
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Figure  C-8.  Earth  (a-7500  km):  Altitude  of  Arc 
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C-9.  Earth  (a=7500  km):  Duration  of  Arc  Ah  <  . 1  km 
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Figure  C-11.  Earth  (a-7500  km):  Latitude  Range  of  Arc  Ah  <  . 1  km 
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Figure  C-12.  Earth  (a-7500  km):  Latitude  Range  of  Arc  Ah  <  1km 
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Figure  C-13.  Mars  (a“3600  km):  Mid-Latitude  of  Arc 


Figure  C-14.  Mars  (a-3600  km):  Altitude  of  Arc 
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Figure  C-16.  Mars  (a-3600  km):  Duration  of  Arc  Ah  <  1km 
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Figure  C-17.  Mars  (a=3600  km):  Latitude  Range  of  Arc  Ah  <  .1  km 
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Figure  C-18.  Mars  (a-3600  km);  Latitude  Range  of  Arc  Ah  <  1km 
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Figure  C-19.  Mars  (a=4000  km);  Mid-Latitude  of  Arc 
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Figure  C-22.  Mars  (a-4000  km):  Duration  of  Arc  Ah  <  1km 
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Figure  C-23.  Mars  (a=4000  km):  Latitude  Range  of  Arc  Ah  <  . 1  km 


fUrs  ( »='t000l<n) :  La-li  tude  Range  wi  th  Consiani  RI  ii  ii  de  ( del  ia  K  <  1 .  kn) 


flKajMENT  OF  PERIflPSIS  (deg) 


Figure  C-24.  Mars  (a-4000  km):  Latitude  Range  of  Arc  Ah  <  1km 


Figures  D-1  through  D-30  present  the  secular  changes  in  semi-major 

axis  (D-l:D-6),  inclination  (D-7;D-12),  argument  of  periapsis 

(D-13:D-18),  eccentricity  (D-19:D-24),  and  longitude  of  the  ascending 

node  (D-25;D-30).  For  each  orbital  element,  the  first  two  figures 

present  results  for  Earth  (a  =  6700  km)  with  ballistic  coefficients  of 

.01  X  10  ^  km^/kg  and  .04  x  10  ^  km^/kg.  The  third  and  fourth  figures 

present  results  for  Mars  (a  =  3600  km)  with  ballistic  coefficients  of 
-  6  2  -  6  2 

.01  X  10  km  /kg  and  .04  x  10  km  /kg.  The  final  two  figures  present 

results  for  Earth  (a  =  7500  km)  and  Mars  (a  =  4500  km)  with  a  ballistic 

”6  2 

coefficient  of  .04  x  10  km  /kg. 
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D-2.  Earth(a-6700kiD,  /9-.04e-6  km  /kg):  Semi- major  Axis  Change 
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Figure  D-3.  Mars(a-3600kin,  Ole-6  km  /kg):  Semi-major  Axis  Change 
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Ear-Ui:  a  •  7506  kn  ba'ta  •  .0H*-6  kn'‘2Ag 


Figure  D-6.  Mars(a-4500kin,  ^-.04e-6  km  /kg):  Semi-major  Axis  Change 
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Figure  D-7.  Earth(a-6700km,  Ole-6  km  /kg):  Inclination  Change 


Figure  D-8.  Earth(a-6700kin,  j0-.O4e-6  km^/kg) :  Inclination  Change 
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Figure  D-9.  Mars(a-3600km,  ^-.Ole-6  km  /kg):  Inclination  Change 
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Figure  D-10. 
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Mars(a-3600km,  /9-.04e-6  km  /kg):  Inclination  Change 
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Figure  D-12.  Mars(a-4500km,  j9-.04e-6  km^/kg);  Inclination  Change 
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Figure  D-17.  Earth(a=7500km,  0~.Ohe-6  km  /kg):  Arg  of  Periapsis  Change 
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Figure  D-18.  Mars(a-4500kn),  ^-.04e-6  km^/kg):  Arg  of  Periapsis  Change 
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3-19.  Earth(a-6700km,  Ole-6  km  /kg);  Eccentricity  Change 
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Figure  D-20.  Ear th(a-6 700km,  p-.OUe-b  km  /kg):  Eccentricity  Change 
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Figure  D-21.  Mars(a=3600km,  ;9-. Ole-6  km  Ag) :  Eccentricity  Change 
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Figure  D-22.  Mars(a-3600km,  ^-.04e-6  km  Ag) :  Eccentricity  Change 
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Figure  D-24.  Mars(a-4500kin,  p-.QUe-6  km  /kg):  Eccentricity  Change 
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Figure  D-26.  Earth (a-6 700km,  ^-.04e-6  km  /kg);  Long  of  Asc  Node  Change 
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Figure  D-27.  Mars(a-3600kni,  Ole-6  km  /kg):  Long  of  Asc  Node  Change 
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Figure  D-28.  Mars(a-3600km,  ^-.04e-6  km  Ag):  Long  of  Asc  Node  Change 
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Figure  D-29.  Earth (a-7 500km,  ^-.04e-6  km  /kg):  Long  of  Asc  Node  Change 
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Figure  D-30.  Mars(a-4500km,  ^-.04e-6  km  /kg):  Long  of  Asc  Node  Change 
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Figure  E-5.  Earth(a=7500km,  /3”.04e-6  km  /kg):  Velocity  Impulse 
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Figure  E-6.  Mars(a-4500km,  ^-,04e-6  km^/kg) :  Velocity  Impulse 
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